﻿using System;
using System.Collections.Generic;
using System.Text;
using System.Windows.Forms;
namespace GeoFly
{
    /// <summary>
    /// Horton超渗产流计算类
    /// </summary>
    public class HortonInfil
    {

        public HortonInfil()
        {

        }

        /// <summary>
        /// 设置参数
        /// </summary>
        /// <param name="currow">行标</param>
        /// <param name="curcol">列标</param>
        /// <param name="dSoilW">土壤水分</param>
        /// <param name="dErr">可能是误差</param>
        public void SetGridPara(GridLayerPara g_GridLayerPara, int row, int col, double SP_Sw)
        {
            int SoilID=(int)g_GridLayerPara.g_SoilLayer.Values[row,col];
            m_currow = row;
            m_curcol = col;
            m_dPreSoilW = SP_Sw;
            

            this.m_dK = g_GridLayerPara.SoilTypes[SoilID].Horton_K;
            this.m_dF0 = g_GridLayerPara.SoilTypes[SoilID].SP_Init_F0;
            this.m_dFc = g_GridLayerPara.SoilTypes[SoilID].SP_Stable_Fc;
        }
        /// <summary>
        /// 计算霍顿超渗产流
        /// </summary>
        public double HortonExcessRunoff()
        {
            double dt0 = this.m_dPreSoilW / this.m_dF0;
            //dt0时的土壤含水量
            double dtmpsw = DTempSoilW(dt0);
            //土壤含水量的变化量
            double dthet = Math.Abs(m_dPreSoilW - dtmpsw);
            double dt = 0;
            int epochs = 0;
            while (dthet > 1e-6)
            {
                if (epochs > 1000)
                {
                    throw new Exception("霍顿下渗率计算迭代不收敛。");
                }
                m_dFt = m_dF0 - m_dK * (dtmpsw - m_dFc * dt0);
                dt = dthet / m_dFt;
                dt0 +=dt;
                dtmpsw = DTempSoilW(dt0);
                dthet = Math.Abs(m_dPreSoilW - dtmpsw);
                epochs++;
            } 
            return this.m_dFt;
        }

        /// <summary>
        /// 求时段dt下的土壤含水量变化
        /// </summary>
        /// <param name="dt"></param>
        /// <returns></returns>
        public double DTempSoilW(double dt)
        {
            double dret = 0;
            dret = m_dFc * dt + (1 - Math.Exp(-1 * m_dK * dt)) * (m_dF0 - m_dFc) / m_dK;
            return dret;
        }
        /// <summary>
        /// 土壤水实际下渗量
        /// </summary>
        public double m_dFt;
        /// <summary>
        /// 初始土壤含水量
        /// </summary>	
        public double m_dPreSoilW;


        public double m_dERR;
        /// <summary>
        /// 下渗曲线系数
        /// </summary>
        public double m_dK;
        /// <summary>
        /// 稳定下渗率
        /// </summary>
        public double m_dFc;
        /// <summary>
        /// 起始下渗率
        /// </summary>
        public double m_dF0;
        public int m_currow;
        public int m_curcol;

        
    };

}